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Abstract 



Numerical calculations have been made on the spike-train response of a pair of Hodgkin- 
Huxley (HH) neurons coupled by synapses and axons with time delay. The recurrent 
excitatory-excitatory, inhibitory-inhibitory, excitatory-inhibitory, and inhibitory-excitatory 
couplings are adopted. The coupled, excitable HH neurons are assumed to receive the two 
kinds of spike-train inputs: the transient input consisting of M impulses for the finite dura- 
tion (M: integer) and the sequential input with the constant interspike interval (ISI). The 
HH neurons in all the kinds of couplings are found to play a role of memory storage with on- 
off switching. When the coupling strength and the time delay are changed, the distribution 
of the output ISI To shows bifurcation (multifurcation), metastability and chaotic behavior. 
The calculation of the time correlation shows that the synchronization between the two HH 
neurons is well preserved even when the distribution of their Tq is chaotic. The correlation 
dimension of the cycles of is shown to depend not only on the model parameters but also 
on the type of input ISIs. 
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I. INTRODUCTION 



Neurons communicate by producing sequence of action potentials or spikes. It has been 
widely believed that information is encoded in the average rate of firings, the number of 
action potentials over some suitable intervals. This firing rate hypothesis was first proposed 
by Andrian from a study of frog, in which the firing rate monotonically increases with 
an increase of the stimulus strength. By applying the firing rate hypothesis, the properties 
of many types of neurons in brain have been investigated and the theoretical models have 
been developed 0. 

When all action potentials are taken to be identical and only the times of firing of a 
given neuron are considered, we obtain a discrete series of times, {tn}, which is expected to 
contain the information. In the rate coding, only the average of the rate of the interspike 
interval (ISI) is taken into account, and then some or most of this information is neglected. 

In recent years, the alternative temporal coding, in which detailed spike timing is taken 
to play an important role, is supported by experiments in a variety of biological systems: 
sonar processing of bats [Q, sound localization of owls [Q, electrosensation in electric fish [Q, 
visual processing of cats ||^ monkeys ^ and human It is now primarily important to 
understand what kind of code is employed in biological systems: rate code, temporal code 
or others [|l^ ||TI |. 



Neural functions are performed in the activity of neurons. Since the Hodgkin-Huxley 
(HH) model was proposed to account for the squid giant axon [0], its property has been 
intensively investigated. Its responses to an applied dc ||T^ and sinusoidal currents 
p!8| |1I9| have been studied. The HH-type models have been widely employed for a study 
on activities of transducer neurons such as motor and relay neurons, which transform 
amplitude-modulated inputs to spike-train outputs. Regarding the single HH neuron as 
a data-processing neuron, the present author [^] (referred to as I hereafter) recently in- 
vestigates its response to the spike-train inputs whose ISIs are modulated by deterministic, 
semi-deterministic (chaotic) and stochastic signals. 

Several investigations have been reported on the property of a pair of the HH neurons 
2T| - ||27|| . In the network of two HH oscillators coupled by excitatory couplings without 



time delay, the unit fires periodically in the synchronized state. It is, however, not the case 
when the excitatory couplings have some time delay, for which the antiphase state becomes 
more stable than the synchronized state |^]. Rather, inhibitory couplings with substantial 



time delay lead to the in-phase synchronized states in the coupled HH oscillators [^. The 



similar conclusion is obtained also in the coupled integrate-and-fire (IF) oscillators E^ 
P8|| - [ pl[] . The phase diagrams for the synchronized state and various cluster states in the 
coupled HH oscillators are obtained as functions of the synapse strength and the time delay 

El- 



In recent years, much attention has been paid to the delayed-feedback systems described 
by the delay-differential equation (DDE) Their property has been investigated 

with the use of various functional forms for the delay-feedback term in DDE. The exposed 
properties include the odd-harmonic solutions |Q , the bifurcation (multifurcation) lead- 
ing to chaos ||3^- [0, the multistability [Q, and the chaotic itinerancy ||3^. Among 
them the multistability is intrigue because it may be one of conceivable mechanisms for 
memory storage in biological neural networks. It has been shown by Ikeda and Matsumoto 
that when the delay time is larger than the response time in the delayed feedback system. 
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information may be stored in temporal patterns. Actually, Foss, Longtin, Mensour and Mil- 



ton [25 1 demonstrate this ability in the coupled HH (and IF) neurons with the time-delayed 
feedback. 

Recurrent loops involving two or more neurons with excitatory and/or inhibitory 



synapses are found in biological systems such as hippocampus [Q, neo-cortex |^ and 
thalamus [^. It is important to make a detailed study on the coupled HH neurons, which 



is the simplest but meaningful network unit. The purpose of the present paper is to inves- 
tigate the spike-train responses of the coupled, excitable HH neurons, which is treated as 
data-processing neurons as in our previous study |^|. We consider the four kinds of recur- 



rent couplings between a pair of HH neurons: the excitatory-excitatory (E-E), inhibitory- 
inhibitory (I-I), excitatory- inhibitory (E-I), and inhibitory-excitatory (I-E) couplings. We 
apply, to the coupled HH neurons, the two kinds of spike-train inputs: one is the transient 
inputs consisting of clustered M impulses for the finite duration (M: integer) and the other 
is the sequential spike train with the constant ISI. The transient and asymptotic property 
of their response will be investigated. 

Our paper is organized as follows: In the next Sec. II, we describe a simple neural 
system consisting of neurons, axons, synapses and dendrites, which are adopted for our 
numerical calculation. We present the calculated results in Sec. HI: the response of the 
coupled HH neurons to the transient, clustered impulses is discussed in Sec. IIIA and that 
to the sequential spike-train input in Sec. IIIB. The dependences of the distribution of the 
output ISIs on the coupling strength and the time delay are studied. The final Sec. IV is 
devoted to conclusion and discussion. 

II. ADOPTED MODEL 

We adopt a simple neural system consisting of a pair of neurons which is numbered 1 and 
2. The neurons which are described by the HH model with identical parameters, are coupled 
with the time delay of Tjk {j, k = 1,2) for an impulse propagating from the neuron k to the 
neuron j. This delay time is the sum of conduction times through the axon and dendrite. 
It has been reported that real biological synapses exhibit temporal dynamics of depression 



or potentiation during neuronal computation |^2|. We, however, treat the synapse as a 



static unit for a simplification of our calculation. The synapse with the coupling strength 
Cjk is excitatory or inhibitory, and it is assumed to be described by the alpha function 
[Eq.(7)]. 

Dynamics of the membrane potential Vj of the coupled HH neuron j (=1, 2) is described 
by the non-linear DDEs given by 

CdV,{t)/dt = -ir(yj,m„ h„n,) + If + If\mt - r,,)}), (1) 

where (7 = 1 /iF/cm^ is the capacity of the membrane. The first term of Eq.(l) expresses 
the ion current given by 

lf{V,, m„ h,, n,) = gn.m]h,{V, - V^^) + qku^V, - Vk) + gi.{V, - Vi^. (2) 

Here the maximum values of conductivities of Na and K channels and leakage are ^^Na = 
120 mS/cm^, = 36 mS/cm^ and gi, = 0.3 mS/cm^, respectively; the respective reversal 



3 



potentials are Vj^^ = 50 mV, Vk = —77 mV and Vl = —54.5 mV. The gating variables of 
Na and K channels, m^, hj and described by 



drrij/dt = -{amj + bmj) rrij + amj, (3) 
dhj/ dt = -{tthj + bhj) hj + ahj, (4) 

drij/dt = -{ttnj + bnj) rij + a„j. (5) 

The coefficients of amj and hmj etc. are expressed in terms of Vj (their explicit expressions 
having been given in Refs. |^ [^) and then the variables Vj, rrij, hj and rij are coupled. 
The second term in Eq.(l) denotes the external input currents given by 

/f* = /s, + ^5,i5]a(t-t.„), (6) 

n 

with the alpha function a{t) given by 

ait) = (t/r^) e-*/- e(t), (7) 

The ffist term (Isj) in Eq.(6) is the dc current which determines whether the neuron is 
excitable or periodically oscillating. Its second term expresses the postsynaptic current 
which is induced by the presynaptic spike-train input applied to the neuron 1, given by 

Ui{t)=V,Y.^it-t^n). (8) 

n 

In Eqs.(6)-(8), Q{t) = 1 for x > and for x < 0; = gs{Vs, - K), fi-s and V; stand for 
the conductance and reversal potential, respectively, of the synapse; Tg is the time constant 
relevant to the synapse conduction, which is assumed to be Tg = 2 msec; t^ is the n-th ffiing 
time of the spike-train inputs defined recurrently by 

^in+l i'm ~l~ -^in('^in); (9) 

where the input ISI Ti„ is generally a function of tin- For the constant input ISI of Tin = Ti, 
tin is given by tin = nTi for an arbitrary integer n. 

When the membrane potential of the j-th neuron Vj{t) oscillates, it yields the spike-train 
output, which may be expressed by 

Uojit)=V,Y.^i^-^ojm), (10) 
m 

in a similar form to Eq.(8), tojm being the m-th ffiing time when Vj{t) crosses Vz = 0.0 mV 
from below. The output ISI is given by 

Tojm ^oj'm+l tojm- (H) 

The third term in Eq.(l) which expresses the interaction between the two neurons, is 
assumed to be given by 
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lf\{Vk{t - r,fc)}) = E E (^Jk ait - - t^km)- (12) 

The positive and negative Cjk denote the excitatory and inhibitory couphngs, respectively. 
Depending on the signs of Cjk, we consider the following four types of couplings (Fig.l): 

(a) C21 > and C12 > (E-E coupling), 

(b) C21 < and C12 < (I-I coupling), 

(c) C21 > and C12 < (E-I coupling), 

(d) C21 < and Cu > (I-E coupling). 

In order to reduce the number of parameters, we assume that | C21 | = | C12 \= c and 

As for the functional form of the coupling term of /™*({\4(t — t^^)}), Foss, Longtin, 
Mensour and Milton |^ adopt a simpler form given by 

lf\{Vk{t - r,fc)}) = f'jk Vk{t - T,k), (13) 

taking no account of the synapse, where fijk is the coefficient of the synaptic coupling. 
They discuss the memory storage of the pattern in output spike trains, injecting the input 
information by the initial function, V{t) for t G [—Td,0), whereas in our calculation input 
information is given by JJ^* [Eq.(6)]. 

Differential equations given by Eqs.(l)-(5) including the external current and the coupling 
given by Eqs.(6)-(12) are solved by the forth-order Runge-Kutta method. The calculation 
for each set of parameters is performed for 2 sec (200,000 steps) otherwise noticed by the 
integration time step of 0.01 msec with double precision. The initial conditions for the 
variables are given by 

Vj{t) = -65 mV,mj{t) = 0.0526, /ij(t) = 0.600, ^^(t) = 0.313, for j = 1,2 att = 0, (14) 

which are the rest-state solution of a single HH neuron [cjk = 0). The initial function for 
Vj(t), whose setting is indispensable for the delay-differential equation, is given by 

Vj{t) = -65 mV for j = 1, 2 at t G [-r^, 0). (15) 

For an analysis of asymptotic solutions, we discard results of initial 200 msec (20,000 steps). 



III. CALCULATED RESULTS 

In the present study, we consider only the excitable HH neurons by setting Isj = and 
As = 40/iA/cm^. Our model has additional three parameters, Ti, and c. We treat them as 
free parameters to be changed for the four kinds of E-E, I-I, E-I, and I-E couplings because 
the values of ISI and the time delay observed in biological systems distribute in a fairly wide 
range |l43 |. 



A Transient Spike-Train Inputs 
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Let us first investigate the response to the transient, clustered spike-train inputs consist- 
ing of M impulses. Figure 2(a) shows the example of the time courses of input (fA), output 
pulses {Uoj), the total postsynaptic current (Jj = /|^* + and the membrane potential 
{Vj) with M = 3, Ti = 20, = 10 msec and c = 1.0 for the E-E coupling (c > 0). The 
first external pulse applied at t = yields the firing of the neuron 1 after the intrinsic delay 
of Til ~ 2 msec. The emitted impulse propagates the axon and reaches the synapse of the 
neuron 2 after T21 = 10 msec. After a more delay of an intrinsic Ti2 ~ 2 msec, the neuron 
2 makes the firing which yields the input current to the neuron 1 after a delay of T12 = 10 
msec. The input pulses trigger the continuous oscillation in the coupled HH neurons with 
the output ISI of To = 24.1 msec. Filled and open circles in Fig. 2(b) denote the time 
dependence of the output ISI of the neuron 1 and 2, respectively. We note that Tqi and T02 
start from the values of 20.00 and 19.96 msec, respectively, and soon become the fixed value 
of 24.10 msec. 

Figure 3(a) shows the result for the I-I coupling with the parameters same as in Fig. 
2(a). After T21 = 10 msec of the firing of the neuron 1, the inhibitory input arrives at 
the neuron 2, which fires by the inhibitory rebound process with an additional delay of 
Til ~ 15 msec beside a delay of T12 = 20 msec of the axon. Subsequently the neuron 1 
induced by the recurrent impulse coming from the neuron 2, fires after an total time delay 
of T21 + ~ 25 msec due to delays of the axon and the inhibitory rebound. Figure 3(b) 
shows that the output ISIs of the neuron 1 and 2 start from the values of 20.00 and 21.04 
msec, respectively, and asymptotically approach the oscillating two values of 24.33 and 24.45 
msec. 

In the cases shown in Figs. 2 and 3, the delay of ra = 10 msec is shorter than the 
duration of input pulses (40 msec). On the contrary, when the delay time becomes larger 
than the duration time of clustered inputs, the situation is changed. Figures 4 shows the 
time courses of the membrane potentials of the neuron 1 with the four couplings for M = 3, 
Ti = 20, Td = 50 msec and c = 1.0 (hereafter we consider only the output of the neuron 1). 
We note that the triggered oscillation continues, keeping the original three-impulse form in 
all the couplings. Then the coupled HH neurons play a role of memory storage |3^ . 



The period of the induced oscillation Tp shown in Fig. 4 is essentially determined by the 
total delay time of the feedback loop Tfb given by 

~ Tfb = 2rd + Til + Ti2 (16) 

where Ty denotes the intrinsic time delay of the neutron j. It is noted that Xy depends 
on whether an input is excitatory or inhibitory: is about 2-3 msec for the excitatory 
input while it is about 14-15 msec for the inhibitory input. Our simulation shown in Fig. 
4 yields Tp™ = 105, T^^ = 129, Tp^^ = T^^ = 117 msec for the E-E, I-I, E-I and I-E 
couplings, respectively, which follows Eq.(16) with ti = 2.5 and 14.5 msec for the excitatory 
and inhibitory synapses, respectively (r^ = 50 msec). 

I. The coupling-strength dependence 

Figure 5 shows the coupling-strength dependence of the period and the output ISIs of 
the asymptotic solution in the cases of the four couplings for M = 3, Ti = 20 and Td = 50 



msec 1^^. The period of the oscillation triggered by input pulses has little dependence on 
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c, as far as c is larger than the critical value, Ccr, above which the oscillation continues. In 
the case of M = 3, Tj = 20 and Td = 50 msec, for example, we get Ccr = 0.20 for the E-E 
coupling and 0.42 for the I-I, E-I and I-E couplings: note that at 0.44 < c < 0.76 in the 
latter cases, we get the two-impulse output against the three- impulse input. This is because 



the firing by the inhibitory rebound requires a considerable synaptic magnitude ^0 

2. The time-delay dependence 

Next we study how the output ISIs are determined. When the feedback time Tfb is larger 
than the duration of clustered impulses {i.e. Tfb = 2rd + + > (M — 1) T{), we get two 
values of Tq given by 

^ y.^ j.i2) ^ _ - 1) Ti = 2rd + ra + r,2 - (M - 1) T,. (17) 

On the other hand, when the feedback time is shorter than input-pulse duration (2rd + + 
< (M-l)Ti), we get 

T« = Ti e(M - 3), Ti^) =1 £ Tfb - m Ti 1 = 1 £ (2rd + + t,^) ~mT,\, (18) 

where integers i and m satisfy 1 < £ < [{M — l)Ti/Tfb] + 1 and < m < M — 1, [ ■ ] is the 
Gauss sign and T^^^ is vanishing for M < 2. 

We have obtained the ra dependence of the output ISI from the asymptotic solution in 
our simulations. Results for clustered inputs of M = 3 are shown in Fig. 6 (a)-6(d), where 
filled circles denote the values of Tq for a given [Q. Figure 6(a) shows the result for 
inputs with Ti = 20 msec in the E-E couplings. We try to analyze the calculated Tq with the 
use of Eqs.(17) and (18). Three dashed lines in Fig. 6(a) are given by 2rd -|- 5, 2rd — 15 and 
2rd — 35, which are obtained by Eqs.(17) and (18) with th = t\2 = 2.5 and Tj = 20 msec. 
Dashed lines in Fig. 6(b) are given by 2rd + 29, 2rd + 9 and 2rd — 11, which are obtained 
by Eqs.(17) and (18) with th = = 14.5 and T[ = 20 msec. Similar analysis is made for 
the 15-1 and I-E couplings, whose results are expressed by dashed lines in Figs. 6(c) and 
6(d). Agreement of the dashed lines with filled circles in Fig. 6(a) and 6(b) are fairly good. 
On the contrary, the agreement between them is not satisfactory in Fig. 6(c) and 6(d) for 
Ti = 50 msec. The deviation of filled circles from the dashed lines is due to the effects of 
non-linearly of the HH neurons. 

The time courses of the membrane potentials for M = 1 — 5 are plotted in Fig. 7 for 
Ti = 20, Tci = 50 msec and c = 1.0. In all the couplings, the form of input pulse is preserved 
and the coupled neurons play a role of memory storage. 

It is possible to control the switching of the oscillations by input pulses. Figure 8 demon- 
strates such a switching in the E-I coupled HH neurons. The oscillation is triggered by the 
external input Ui at t = 0. When we apply the external excitatory input at t = 543 msec, 
it cancels out the recurrent inhibitory input to the neuron 1 and suppress the firing of the 
neuron 1. Then the oscillation is switched off by the control input pulse at t = 543 msec. At 
t = 700 mec, the oscillation is again triggered by input pulse. In this example, the control 
pulse is applied to the neuron 1. We can alternatively switch off the oscillation by applying 
the inhibitory (excitatory) input to the neuron 2 in the E-I (I-E or I-I) coupling. These show 
the feasibility of the on-off switching of the oscillation in the simplest neuronal unit. 
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B. Sequential Spike-Train Inputs 



Next we discuss the response to the sequential spike-train. Our calculations in 1 show 
that when an isolated HH neuron (c = 0) receives the sequential inputs with the constant ISI 
of Ti, it behaves as a low-pass filter: it emits the spike train with Tq > 10 msec for T[ < 12 
msec while for Tj > 12 msec its output ISI is given by Tq = (see Fig. 7 of Ref. ||2^). This 
response may be modified when the coupling is introduced to a pair of HH neurons. Figures 
9(a) and 9(b) show the time courses of input (f/i), output (Uoj), the total postsynaptic 
current {Ij = JJ^* + Jj"^*) and the membrane potential [Vj) for Ti = 20 msec, = 10 msec 
and c = 1.0 for the E-E and I-I couplings, respectively. They should be compared with Figs. 
2(a) and 3(a) for the transient spike-train inputs with the same values of and r^. The 
output ISI in Fig. 2 is 24.1 msec while that in Fig. 9(b) is 20 msec which is the entrained 
value with input ISI. Although To in Fig. 3 are two values of 24.3 and 24.5 msec, we get, in 
Fig. 9(b), ten values of 20.6, 21.2, 37.3, 19.2, 22.1, 20.6, 21.7, 35.7, 19.9 and 21.8 msec. 

I. The coupling-strength dependence 

The response behavior of the coupled neurons strongly depends on the parameters of c, 
Td and Ti. Figure 10(a)-10(d) show the examples of the c dependence of the distribution 
of To of the E-E, I-I, E-I and I-E couplings for various parameters We note that, as 
increasing the c value, the distribution of the output ISIs show the bi(multi)furcation, as 



commonly observed in systems with the delayed feedback ||33|. In order to investigate the 
phenomenon in more detail, we show, in Fig. 11(a), the enlarged plot for the range of 
0.6 < c < 1.2 sandwitched by the dotted, vertical lines in Fig. 10(a). Figure 11(b) is the 
enlarged plot of Fig. 10(d) for the narrow range of 0.3 < c < 1.0. These figures clearly show 
the bi(multi)furcation with many windows. 

A cycle whose output ISIs almost continuously distribute, is expected to be chaotic 
although in the strict sense, the distribution of ToS never becomes continuous because they 
are quantized by the integration time step of 0.01 msec. Among many candidates of chaos- 
like behavior in Figs. 10 and 11, we pay our attention to the result of c = 0.95 in Fig. 11(a), 
for which the Lorentz plot (return map) of its To is shown in Fig. 12(a) (calculations are 
performed for 20 sec of two million steps). The output ISIs seem to distribute on the folded 
ring. When these points are connected by lines in the chronological order, the inside of the 
ring is nearly filled by them. In order to examine the property of this cycle, we calculate 
the correlation dimension u given by 



, log C(e) , , 

lim f ^ \ (19 
log e ^ ^ 



with 

N 



C(e) = J2 \Xrr^-X^ |), (20) 



m.n=l 



■^m {Torri) ; -^m+fc— l); (^-l) 
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where C(e) is the correlation integral, X^n is the fc-dimensional vector generated by Tom 
N the size of data, and 9(-) the Heaviside function. Figure 12(b) shows the logC(e)-loge 
plot for various k calculated for the cycle shown in Fig. 12(a) with ~ 1200. We note 
that C(e) behaves as C(e) cx with the correlation dimension of u = 0.94 ± 0.02 for small 
e (0.01 = e~^'^ < e < e°). The non-integral u implies that this cycle may be chaos, the 
related discussion being given in Sec. IV. 

2. The time-delay dependence 

A simple generalization of the discussion given in the previous Sec. IIIA2 yields that the 
output ISIs for the sequential spike-train input with the constant T; are given by 

To =1 £ Tfb - m Ti | = K (2rd + Til + ^2) - m T; |, (22) 

where i, m are integers. In the above analysis, we do not take into account the absolute 
refractory period and the merging of recurrent inputs with external inputs at synapse, which 
eliminate some of the candidates of Tq obtained from Eq.(22) for many choices of a pair of 
i and m, as will seen shortly. 

Figure 13(a)-13(d) show the dependence of the distribution of Tq in the asymptotic 



solution for the sequential inputs |4J]. We try to analyze the obtained distribution of Tq 
with the use of Eq.(22) by adopting Tij=2.5 and 14.5 msec for the excitatory and inhibitory 
couplings, respectively, as was made in Sec. IIIA2. Dashed curves in Fig. 13(c)-13(d) denote 
the lines generated by Eq.(22) for these values of Ty and Tj = 50 msec with a choice of a 
pair of integers {i, m) shown beside the lines. Although many lines can be drawn for various 
choices of the integers, we plot only results with £ = 1 in order to avoid a disfigurement by 
them. Some of filled circles representing the distribution of Tq in Figs. 13(c)-13(d) may be 
explained by dashed lines. Most of filled circles, however, cannot be well understood by this 
analysis due to the effects not included in our simple analysis. Furthermore our analysis with 
the use of Eq.(22) does not work at all for the results shown in Fig. 13(a) and 13(b), which 
plot the distribution of Tq in the E-E and I-I couplings for inputs with T[ = 20 msec. We 
note that the result shown in Fig. 13(a) (13(b)) for the sequential inputs is quite different 
from that shown in Fig. 6(a) (6(b)) for the transient inputs although both the calculations 
adopt the same parameters of Ti and c. 

In order to see more the detailed structure of the bi(multi)fucation, we show, in Fig. 
14(a), the enlarged plot for the range of 21 < < 26 msec between the dotted, vertical 
lines in Fig. 13(a). Figure 14(b) is the enlarged plot of Fig. 13(c) for the narrow range of 
33 < Td < 38 msec. They clearly show the bi(multi) furcation as changing rj. 



VI. CONCLUSION AND DISCUSSION 

We have performed numerical investigation on the spike-train responses of a pair of HH 
neurons with four kinds of the recurrent E-E, I-I, E-I and I-E couplings. The response of 
the coupled, excitable HH neurons to the transient, clustered impulses and to the sequential 
inputs shows a rich of variety. The former input triggers the synchronous oscillation in the 
coupled HH neurons, which may play a role of memory storage. The response to the latter 
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sequential input strongly depends on the coupling strength and the time delay, yielding 
bi(multi)furcation, multistabilty and chaotic behavior. 

We have applied the two types of inputs of the transient and sequential spike-train 
impulses. On the theoretical point of view, the latter is taken as the limit of M — oo of 
the former. Figure 15(a) and 15(b) show the T; dependence of the distribution of To for 
the transient and sequential inputs, respectively, with = 50 msec and c = 1.0 in the E-E 
coupling. Dashed curves in Fig. 15(a) are given by the equations shown beside the lines 
which are obtained from Eqs.(17) and (18), as mentioned in Sec. IIIA2. On the other hand, 
dashed lines in Fig. 15(b) are obtained with the use of Eq.(22) for a pair of integers (£, m) 
shown in the brackets. Note that dashed lines given by T[, iT; =p 105 and ±2T[ =p 105 in Fig. 
15(a) are identical with those in Fig. 15(b) with (£, m) = (0,l), (1,1) and (1,2), respectively. 
It is apparent that the distribution of Tq in Fig. 15(a) is not the same as that in Fig. 15(b), 
but they are partly similar. For example, we obtain for T[ = 20 msec, Tq = 19.6, 20.0 and 
64.54 msec in Fig. 15(a) whereas only 20.0 msec in Fig. 15(b). In order to understand 
this difference, we plot, in Fig. 16, the time dependence of Tq for this set of parameters by 
changing the M value. For M = 3, To oscillates with the values of 19.6, 20.0 and 64.54 msec, 
as mentioned above. The calculated To for M=4 are 18.7, 19.8, 20.0, and 45.6 msec, and 
those for M = 5 are 19.9, 20.0 and 24.2 msec. For M = 10, To remains 20 msec until t ~ 200 
msec, after which Tq oscillates with the values of 19.9, 20.0 and 24.2 msec. In the limit of 
M — ^ oo corresponding to the sequential inputs, the state with Tq = 20 msec continues from 
t = to oo. Thus as increasing M, the time region of To = 20 msec is increased. 

Figure 17 shows the similar plot of the time dependence of the distribution of To for 
various M with Tj = 20, ra = 13.75 and c = 0.95, for which the sequential input leads to 
the chaotic behavior, as was discussed in Sec. IIIBl (see Fig.l2(a)). In the case of M = 3, 
we get the oscillation in Tq which asymptotically approaches the value of 15.97 msec. In the 
case of M = 10 (50), the chaotic behavior is realized at < t ~< 180 {0 < t ~< 980) msec 
during the application of inputs. After inputs are switched off, the output ISI gradually 
approach the fixed value of 15.97 msec. In the limit of M ^ oo, the chaotic oscillation 
eternally continues. 

We have shown in Sec. IIIBl that the cycle of the output ISIs shown in Fig. 12(a) may 
be chaos because its correlation dimension of z/ ~ 0.94 is derived from the log C(e)-log e plot 
in Fig. 12(b). This is not surprising because the response of single HH neurons to some 
kinds of external inputs may be chaotic |jl8| |T9| In particular, it has been shown in I 



that the response of a single HH neuron may be chaos when the ISI of the spike-train input 
is modulated by the sinusoidal signal [pO[] : 



Ti(t) =do + di sin(27rt/Tp). (23) 

Figure 18(a) shows the Lorentz plot of the output ISIs of the single HH neuron receiving 
sequential inputs with ISIs sinusoidally modulated by Eq.(23) with do = 2di = 20 and 
Tp = 100 msec (see Fig. 7(d) of I). We note that ToS distribute on the deformed ring. From 
the logC(e)-loge plot (not shown) of this cycle, we get its correlation dimension of ~ 1.04. 
We apply this spike-train input to the E-E coupled HH neurons with Td = 10 msec and 
c = 1.0, whose Lorentz plot is shown in Fig. 18(b). Its structure is rather different from 
that shown in Fig. 18(a). Actually the correlation dimension of this cycle for the coupled 
HH neurons is ~ 1.83, which is different from and larger than u ~ 1.04 of the cycle shown 
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in Fig. 18(a) for the single HH neuron. From similar calculations for the coupled HH nurons, 
we obtain the correlation dimensions of z/ ~ 0.95 for Td = 5 msec and c = 1.0, and v ~ 1.03 
for Td = 10 msec and c = 0.5. These results clearly show that the correlation dimension 
of the output ISIs depend not only on the model parameters (c and Td) of the coupled HH 
neurons but also on the correlation dimension of input ISIs {v\ = for the constant ISI and 
z/j = 1 for the sinusoidally modulated ISI). We expect that spike-train inputs with larger 
correlation dimension lead to spike-train outputs with larger u. One of the disadvantages of 
the present calculation of the correlation dimension is a lack of the data size of ~ 1200 
with million-step calculations. A more accurate analysis requires a larger size of data and 
then a computer with the larger memory storage. 

Next we discuss the time correlation ri2(r) between the membrane potentials, Vi and 
V2, of the neurons 1 an 2, defined by 



where ta = 1000 and U = 2000 msec are adopted for our calculation. Figure 19(a) shows 
the result for the case of the sequential input to the E-E coupled HH neurons with Ti = 20, 
Td = 10 msec and c = 1.0 (see Fig. 9(a) for the time courses of Vi and V2). In this case we 
obtain the constant Tq = 20 msec as was discussed in Sec. IIIB, and then ri2(T) shown in 
Fig. 19(a) has peaks at r = 12.04 + 20n msec (n: integer) with the period of 20 msec, as 
expected. We are interested in the time correlation for the case when the distribution of Tq 
is chaotic. Results for such cases are shown in Figs. 19(b) and 19(c). We have discussed 
in Sec. IIIBl that the cycle of Tq depicted in Fig. 12 may be chaotic. Figure 19(b) shows 
the result of this case for the E-E coupling with Tj = 20, Td = 13.75 msec and c = 0.95. 
We note that ri2(r) has peaks at r=0.0, 16.07, ~ 3.2, 48.17, 62.71,... msec with the period 
of about 16 msec, which is the sum of and th. More evident peaks are found in Fig. 
19(c) showing also the chaotic case discussed in the preceding paragraph: the E-E coupled 
neurons receiving the sinusoidal inputs given by Eq.(23) with do = 2di = 20, Tp = 100, 
Ti = 20, Td = 10 msec and c = 1.0 [see Fig. 18(b)]. We note peaks in ri2(T) at r=0.0, 
12.72, 25.30, 37.88, 50.43, ... msec with the period of about 12.6 msec. It is interesting 
that the synchronization is well preserved between the coupled HH neurons even when the 
distribution of their output ISIs shows the chaotic behavior. 

A fairly large variability (cv = 0.5 ~ 1.0) has been reported for spike trains of non- 
bursting cortical neurons in VI and MT of monkey [Q. It is possible that when the 
appreciable variability in neuronal signals is taken into account in our calculations, much of 
the fine structures in the c— and Td-dependent distributions of will be washed out. In 
order to study this speculation, we apply the spike-train input with ISI whose distribution 
is given by the gamma distribution defined by pO|| 



where F (r) is the gamma function. The average of input ISI is given by fii = r/s, its 
root-mean-square (RMS) by cXi = y/r/s and its variability by Cvi = (Ji/ fii = 1/v^- Figure 
19 shows the dependence of the mean (/Xq) and RMS values ((Tq) of the output ISIs for 
Cvi = 0.0 (dashed curves) and c^i = 0.43 (solid curves) with Ti = 20 msec and c = 1.0. Note 
that cTq provides us with the measure of the width of the distribution of Tq. The distribution 




(24) 



P(T) = T'-^ e-'^l r(r) 



(25) 
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for Cvi = has a fine structure reflecting the strong ra dependence of To [see Fig. 13(a)]. 
This fine structure is, however, washed out for Cyj = 0.43, as expected. Detailed calculations 
of the response of the coupled HH neurons to stochastic spike-train inputs are now under 
way and will be published elsewhere. 
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FIGURES 



FIG. 1. Adopted pairs of HH neurons with (a) the E-E, (b) I-I, (c) E-I, and (d) I-E couphngs, 
open and filled circles denoting excitatory and inhibitory synapses, respectively. 

FIG. 2. The time dependence of (a) the transient input (Ui), output (Uoj), the total postsy- 
naptic current (Ij) and the membrane potential (Vj), and (b) the output ISI (Toj) with M = 3, 
Ti = 20 and = 10 msec in the E-E couplings. The result of V2 in (a) is shifted downward by 200 
mV and scales for Ui, Uoj and Ij are arbitrary. 

FIG. 3. The time dependence of (a) the transient input (U), output (Uoj), the total postsy- 
naptic current (Ij) and the membrane potential {Vj), and (b) the output ISI (Toj) with M = 3 
and Ti = 20 msec in the I-I couplings. Same as in Fig. 2. 

FIG. 4. The time course of the input (L^i) and membrane potentials {V) for M = 3, Ti = 20 
and Td = 50 msec in the four couplings, results of I-I, E-I and I-E couplings being shifted downward 
by 200, 400 and 600 mV, respectively. 

FIG. 5. The coupling-strength dependence of the oscillation period Tp and the output ISIs Tq 
of the E-E (solid curve, open circles), I-I (dashed curve, filled circles), and E-I and I-E couplings 
(dotted curve, open triangles) for the transient inputs with M = 3, Ti=20 and rd=50 msec. The 
left (right) ordinate is for Tq (Tp). 

FIG. 6. The Td dependence of the distribution of To of (a) E-E (Ti = 20 msec), (b) I-I (Ti = 20 
msec), (c) E-E (Ti = 50 msec), and (d) I-E couplings (Ti=50 msec) for the clustered inputs with 
M = 3 and c = 1. The dashed lines are expressed by the equations written beside the lines (see 
text). 

FIG. 7. The time course of the membrane potentials (Vi) for M = 1 — 5 for Ti = 20 and rd=50 
msec with the E-E coupling. 

FIG. 8. The time course of the external input (Ui), the postsynaptic current (h), and the 
membrane potential (Vi) in the E-I coupling with rd=50 msec, demonstrating the switching of the 
oscillation hy Ui. Scales for Ui and Ii are arbitrary. 

FIG. 9. The time course of input (U), output (Uoj), the total postsynaptic current (Ij), and 
the membrane potential (V^) for sequential input with Ti = 20 and Td = 10 msec in the (a) E-E 
and (b) I-I couplings. 
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FIG. 10. The distribution of To as a function of the coupling strength c of (a) E-E (Ti=20, 
rd=13.75 msec), (b) I-I (Ti=20, rd=20 msec), (c) E-I {Ti=20, rd=15 msec), and (d) I-E couphngs 
(ri=20, rd=20 msec) for the sequential inputs. The enlarged plots of the regions between dotted, 
vertical lines in (a) and (c) are shwon in Fig. 11(a) and (b), respectively. 

FIG. 11. The enlarged plots of the distribution of To for (a) the E-E coupling with T = 20 
and Td = 13.75 msec (see Fig. 10(a)), and for (b) the I-E coupling with T = 20 and Td = 20 msec 
(see Fig. 10(c)). The arrow in (a) denotes the c value for which the Lorentz plot is shown in Figs. 
12(a). 

FIG. 12. (a) The Lorenz plot of To for c = 0.95 with T[ = 20 and Td=13.75 in the E-E coupling, 
the computation being performed for 20 sec (two million steps), (b) The correlation integral C(e) 
of the cylce shown in (a) as a function of e in the log-log plot for various dimensions k, the dashed 
line denoting C(e) oc with the correlation dimension of = 0.94 (see text). 

FIG. 13. The Td dependence of the distribution of To of (a) E-E (Ti = 20 msec), (b) I-I 
(Ti = 20 msec), (c) E-E (Ti = 50 msec), and (d) I-I couphngs (T; = 50 msec) for the sequential 
input with c = 1.0. Dashed lines in (c) and (d) are given by Eq.(22) with a pair of integers {i,m) 
shown beside the line. The enlarged plots of the regions between dotted, vertical lines in (a) and 
(c) are shown in Fig. 14(a) and (b), respectively. 

FIG. 14. The enlarged plot of the distribution of To for (a) Ti = 20 msec (see Fig. 13(a)) and 
(b) Ti = 50 msec (see Fig. 13(c)). 

FIG. 15. The Ti dependence of the distribution of Tq for (a) the transient (M = 3) and (b) 
sequential spike-train input with Td = 50 msec and c = 1.0. Dashed lines are given by Eqs.(17), 
(18) and (22) (see text). 

FIG. 16. The time dependence of To for the clustered impulse inputs with M = 3, 10, 50 and 
GO with Ti = 20, Td = 50 msec and c = 1.0 in the E-E coupling, results of M=3, 10 and 50 being 
shifted upward by 30, 20 and 10 mV, respectively. The arrows denote the time below which the 
inputs are continuously applied. 

FIG. 17. The time dependence of Tq for the clustered impulse inputs with M = 3, 10, 50 and 
GO with Ti = 20, Td = 13.75 msec and c = 0.95 in the E-E coupling. Same as in Fig.l6. 

FIG. 18. The Lorenz plots of To of (a) the single HH neuron and (b) the E-E coupled HH 
neurons (rd = 10 msec and c = 1.0) receiving spike-train inputs whose ISIs are modulated by 
sinusoidal signal given by Eq.(23) with do = 2di = 20 and Tp = 100 msec (see text). 
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FIG. 19. The time correlation ri2(r) between the membrane potentials of the neurons 1 and 
2 for (a) the constant-ISI input with Ti = 20, = 10 msec and c = 1.0, (b) that with Ti = 20. 
Td = 13.75 msec and c = 0.95, and (c) the sinusoidal input given by Eq. (23) with di = 2d2 = 20, 
Tp = 100, Td = 10 msec and c = 1.0. The results of (b) and (c) are shifted downward by 1.0 and 
2.0, respectively (see text). 

FIG. 20. The dependence of the mean (/Xo) and rms (ctq) of output ISIs of the E-E copuled 
HH neurons (c = 1) rccicving sequential inputs of Tj = 20 msec with Cyi = 0.0 (dashed curves) and 
Cyi = 0.43 (solid curves) (see text). 
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